Extreme value distribution for weakly correlated fitnesses in block model 
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We study the limit distribution of the largest fitness for weakly correlated and identically dis- 
tributed random fitnesses. A fitness variable is obtained by taking a linear combination of a fixed 
number of independent random variables drawn from the same parent distribution and two fitnesses 
are correlated if they have at least one common term in the respective sum. We find that for certain 
class of parent distributions, the extreme value distribution for correlated random variables can be 
related either to one of the known limit laws for independent variables or the parent distribution 
itself. For other cases, new limiting distributions appear. The conditions under which these results 
hold are identified. 



I. INTRODUCTION 

Extreme value theory [l[ deals with the smallest or the largest of a set of random variables {x\, xjv} and has 
found numerous applications in physics j2], engineering Q, biology [|| and finance If the x^s are independent 
and identically distributed (i.i.d.) according to a continuous distribution p(x), the probability of the kih maximum is 
given by [l[ 
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where the cumulative distribution q(x) — Jjf dyp(y) and B(m 7 n) is the beta function. A classic result in extreme 

value theory of i.i.d. random variables states that for large N (and k <C N), the distribution Pjy (x) is of the following 
scaling form [l|: 



(2) 
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In the above expression, while the location factor 6jv and the scale factor ajv depend on the details of the parent 
distribution p(x), the scaling function is determined by the large x behavior of p(x). For k > 1, the distribution 

Pn ( x ) can ^ e reia ted to Pff by a transformation and the maximum value distribution itself can be one of the 
following three types [l|: 

1. Frechet distribution if p(x) decays as a power law 

2. Weibull distribution if p{x) is bounded above 

3. Gumbel distribution if p(x) is unbounded above and decays faster than a power law 

Much less is known about the extreme value statistics when the random variables are not independent. It is 
interesting to investigate how the correlations affect the limiting distributions of the three i.i.d. classes mentioned 
above. If the correlations are strong, one may expect the extreme value distribution to be different from the i.i.d. limit 
laws and universal distributions may not even exist For weakly correlated variables, on the other hand, the i.i.d. 

results may still hold. In this article, we study the extreme value statistics of random variables with weak correlations 
and identify the relevant parameters depending on which the extreme value distribution is seen to be either one of the 
suitably rescaled i.i.d. limit distributions or the (rescaled) parent distribution or a completely different distribution 
unrelated to these. 

Extreme value statistics of correlated random variables has been studied for stationary Gaussian process with 
correlations decaying logarithmically or faster M] and specific physical examples such as random energy model with 
correlated random potential H,Q, 1// noise |lfj| |. directed polymer on Cay ley tree fluctuating interfaces [T3| and 
mass transport models [l3T ] . Here we study this problem in the context of biological evolution where the question of 



'Electronic address: jain@jncasr.ac.in 



2 



000101000010 



000101001010 



FIG. 1: Block model: A binary sequence of length L = 12 is divided into two blocks of equal length I = 6. The sequences on 
the left and right differing by single mutation have correlated fitnesses as the first block is common between the two sequences. 



extreme values naturally arises as a large population evolves to maximise its fitness. During the first step in evolution, 
an individual in an initially unfit population may acquire one or more mutations. If the total population is infinitely 
large, the subpopulation with D mutations is dominated by the fittest mutant in this subpopulation [14j . As several 
experiments have shown that the fitnesses are not completely random (l~5l - fT7| . we are led to study the extreme value 
statistics of correlated fitnesses. 

In order to examine the implications of correlations amongst fitnesses, several models such as NK model [l8| , block 
model [ItJ and rough Mt. Fuji-type model [l9| have been proposed. In this article, we will work with the block 
model which has been employed in recent years to address questions pertaining to evolutionary dynamics [20l - [22l |. 
The extreme value statistics in the block model has been studied for both strongly and weakly correlated fitnesses 
for an exponentially decaying parent distribution p(x) 23]. However the question of the universality of the extreme 
value distributions was not addressed. Here we extend the previous study by considering a broad class of parent 
distributions and focus on weakly correlated but identically distributed fitnesses. 

The rest of the article is organised as follows. The block model of correlated fitnesses is defined in Sec. [H] The 
extreme value distributions are obtained when a single mutation is present in the initial fitness in Sec. 1111 Al and in 
the presence of multiple mutations in Sec. 1111 B[ Finally a discussion and summary is given in Sec. IIVI 



II. BLOCK MODEL OF CORRELATED FITNESSES 



Many biomolecules such as proteins, antibodies and enzymes have natural domains or partitions and can be modeled 
as a sequence of length L divided into several blocks [l7j . In the simplest setting, a sequence represented by a binary 
string of zeros and ones is divided into B blocks of equal length I = L/B , 1 < I < L (see Fig. [T]). A block 
configuration is assigned a random fitness value regardless of its position in the sequence. These block fitnesses are 
chosen independently from a common distribution p(f) which is nonzero on the interval [a, b] and zero elsewhere. 
Then the sequence fitness is obtained by averaging over the fitnesses of the blocks in the sequence. 

If two sequences have several blocks in common, their fitnesses will also be similar and hence correlated. An 
attractive feature of the block model is that the fitness correlations can be tuned with the block length £. For I = 1, 
as two distinct sequences can have up to L — 1 blocks in common, sequence fitnesses are maximally correlated while 
for i — L, we obtain the model with maximally uncorrelated fitnesses as no common blocks are possible [l4j]- To 
see how fitness correlations vary with block length, consider a set of sequences carrying one mutation relative to the 
sequence with all zeros. Then as a single mutation in this sequence leaves B — 1 blocks unchanged, the fitness Wj of 
the one mutant neighbor of the initial sequence is given by 



(B-l)f + fj 
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where /o and fj denote the fitness of the block configuration with all zeros and with £ — 1 zeros and a one at the jth 
locus respectively. Using the fact that the block fitnesses are independently distributed, we find that the correlation 
amongst the sequence fitnesses {wj} is given by 



(WiWj) - {Wi)(Wj} 
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where a 2 is the variance of the block fitness distribution. It is clear from the last equation that the correlations decay 
as the block length I increases towards L. In this article, we are interested in weakly correlated fitnesses which are 
obtained when L — > oo, I — > oo with B fixed. 



III. DISTRIBUTION OF THE LARGEST FITNESS 



Starting from an initial unfit sequence composed of identical blocks (say, all zeros) , we are interested in finding the 
distribution of the largest fitness amongst (^) sequences with D mutations relative to the initial sequence. In this 
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article, we focus on the extreme value distribution for nonindependent and identically distributed (ni.i.d.) random 
variables (for some results on nonindependent and nonidentically distributed fitnesses, see [H, [HI). Such ni.i.d. 
sequence fitnesses are obtained when either the number of mutations D — 1 for any B > 1 or the number of blocks 
B = 2 and D is odd. In the first case (discussed in Sec. IIII A[) . the distribution of sequence fitness Wj defined by © 
is given by 

ProbK) = J™ J™ dfodfjpiMf^S (wj - ( J? - 1 )/o+./^ (5) 

which is independent of j and thus the sequence fitnesses with single mutation are identically distributed. In the 
second case when B — 2 and a sequence has D mutations, the sequence fitness is obtained by averaging over the 
fitness of the first block with d mutations and the second block with d 1 = D — d mutations (see Sec. IIII B]) . On writing 
down the sequence fitness distribution similar to (J5]) , it is readily verified that only those sequences with different block 
configuration in first and second block have identically distributed fitness. Since a sequence with same configuration 
in both blocks can occur only for even D, it follows that ni.i.d. fitnesses are obtained when D is odd. In the rest of 
the article, we assume L to be an even integer and consider D < L/2 as the results for D > L/2 can be obtained on 
replacing D by L — D. 



A. Single mutation in the initial sequence 

We first consider the extreme value statistics of the fitness set {wj} defined by For a given /o, the sequence 
fitness v = ((B — l)/o + f)/B is the largest amongst the set {wj} if / = max{/i, ft}. But as the block fitnesses are 

i.i.d. random variables, this event occurs with a probability P> (f)- Then the probability P^\w) that the largest 
fitness in the set {wj} has a value w can be written as 

/■oo /* oo 

In general, the probability P^ (w) that the fcth maximum has a value w is given by 

/>oc />oo 

ifV) = J J df dfp(fo)PtHf)S(w-v) 



pr(w)= / df dfp(fo)pr(f)s(w-v) (6) 

Jo Jo 
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For large £, using ^ we obtain 

The behavior of the distribution function in the above equation can be classified as follows: 
(i) If at diverges with £, it is useful to rewrite © as 

(k), . f°° 1 ~m /(I -B)z + Bw-b e \ 
Pf'{w) ra B / dz P (z) —F {h) ± '—Z i (9) 



— oo 



In the limit £,w — >■ oo, the ratio z/at in the argument of F^ k ' above can be ignored and we can write 

pM{v,)K*Fw( Bw -~ bt ) (10) 
at \ at J 

where we have used that p(z) is normalised to unity. Thus up to a rescaling, the extreme value distribution function 
for correlated variables falls in the same universality class as the i.i.d. ones if the scale factor at diverges. The scaling 
variable in (flU|) indicates that the distribution decays faster for correlated variables (B > 1) as one would intuitively 
expect. An example of a class of block fitness distribution for which at diverges with I is p(f ) — 5f s ~ 1 e~f 5 . In 
this case, the Gumbel scaling function F^(y) = e~ v e~ e v with the location factor be — (\n£)i and the scale factor 
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FIG. 2: Scaled distribution for the first maximum when the block fitness distribution p(f) = e~^/(2y/f) (main) and e - ^ 
(inset) for B = 2. The points refer to exact integration of © for £ = 10 2 (+), 10 5 (x) and 10 10 (o). The solid lines show CEOj in 
main and (|14[1 in inset. The inset also shows the Gumbel distribution for comparison (dotted curve). 




at — (\n C)~~^~ which diverges with I for < 5 < 1 1]. A good agreement is seen between the exact distribution 
and asymptotic result (JTUJ) in Fig. [2] for 5 = 1/2. 

(ii) If the scale factor ag vanishes as £ — ► oo, a change of variables in © for B > 1 gives 

rl k) M « f ( ~ atZ + ~ b '- Bw ) (ii) 

(12) 

where the last expression is obtained in the limit t, w oo. Thus the probability distribution of the fcth maximum 
is a rescaled parent distribution and independent of k if ag goes to zero with increasing £. Figure [3] shows that for 
Gaussian distributed block fitnesses, the asymptotic distribution in (fl~2]) approaches the exact distribution in as I 
increases. 

(iii) If the scale factor de is independent of £ for large £, one may expect new limiting distributions to arise. For 
example, for exponentially distributed block fitnesses (5 = 1), the scale factor is unity and the limiting distribution 
of the fcth maximum for i.i.d. variables is given by [l| 

jf> ( x) « f<» (y) = C ^ { ;: r \- ky) ~\ e :%> y Z°° ^ 

[k — lj! ^ e e , y — > — oo 

where ?/ = x — In £. Using this result in © for large I, we can write 

; (w) ?s — / dfe B-t e 



B-Uo (fc-1)! 
73 (£e- B "')-^ ^ 



B-l (fc-1) 



dfe- ^ 



ff k -: 



where the last expression is obtained in the scaling limit £, w — > oo keeping u = Bw — bx£ finite and T(a,z) is the 
incomplete gamma function 24]. When u — > — oo, an asymptotic expansion for large e~ u shows that P^ k \w) ~ e~ e 
and thus the behavior of the backward tail is similar to that in the i.i.d. case. On the other hand, when u — > oo, 

the incomplete gamma function is well approximated by the complete gamma function T ^fc — g \_^ 0j except for 
k = 1 and B = 2 and thus gives Pj; k \w) ~ e - u /( B - 1 )_ For k = 1 and B = 2, the incomplete gamma function 
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FIG. 3: Scaled distribution for the first maximum (main) and the second maximum (inset) for the block fitness distribution 
p(f) = 2/e" /2 and B = 2. The plot shown in points refer to exact integration of for I = 10 2 (+), 10 5 (x) and 10 10 (o) and 
the lines to fl2|). 



r(0, e ") ~ it — C where C is Euler-Mascheroni constant but the tail of the maximum value distribution decays 
exponentially fast in this case also. A comparison of the distribution P^ 1 \w) with the Gumbel distribution is shown 
in the inset of Fig. [2] 



B. Multiple mutations in the initial sequence 



We now turn to the case when a sequence is composed of two blocks and carries odd D > 1 mutations. For a 
sequence divided into two blocks of equal length, the sequence fitness is obtained by averaging over the fitness of 
the first block with d mutations and the second block with d! = D — d mutations. Although the integer d runs 
from to D, distinct sequence fitnesses are obtained if d < d u = (D — l)/2. As the first block can have Sd = ( d ) 
independent fitnesses and the second block can have Sd' fitnesses, the number of distinct sequence fitnesses is given 
by N D = Ystlo s d s d ' = (1/2) (£) Hi which reduces to £ for D = 1. 

We are interested in the distribution P^(w) of the largest fitness w among such No correlated fitnesses. It is 
convenient to work with the cumulative distribution V$ (w) defined as 

(!) t.„\ - / AJpW 



V^(w)= / dw'P^(w') (15) 
Jo 

which gives the probability that all the No fitnesses are smaller than w. If Q Sd (w) denote the probability that the 
fitness of all the sequences carrying d mutations in the first block and d 1 in the second block is smaller than w, we 
can write 

pnIm = n &■«(«>) ( ie ) 

Our task thus reduces to finding the distribution Q Sd (w) which can be calculated in a manner similar to the single 
mutation case. For a given d, the largest sequence fitness v = (fd + fd')/2 is obtained when both the block fitnesses fd 
and fd' of the first and second block respectively are the largest amongst the set of Sd and Sd< i-i-d. random variables. 
Since the probability Q Sd (w) equals the probability that v < w, we get 



g Sd (w) = I dfdj o duPW(f d )P£)(f d ,)e(w-v) 

i 

dfdP^(f d )Vi 1 J(2w-fd) (17) 

poo 

df d .v£(2w-fd,)P$(f d >) (18) 
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where the distribution V(w) = J™ dw'P(w'). 

In the limit £ — > oo, the distributions Q so (w) and Q Sd {w), d > need to be analysed separately. In the former case, 
since sq = 1 on using ([TJ and ((2]) in (fT8f we have 



J —OO ^SD V ^SO 



(19) 



which is simply the cumulative maximum value distribution for sjj random variables when there is a single mutation 
in the initial sequence and therefore exhibits a behavior similar to that already discussed in Sec. IIII A"l For d > 0, as 
both Sd, Sd' are large for £ ^> 1, using (jSJ in (fT7]) we get 

a , w » r * j_f <» ( lM ( ^zizK,) (20) 



* f^-^-^j (21) 

where J : {w) — f dw'F(w') and the last expression is obtained if d Sd /d Sd , — > as £ ~ > oo. The distribution G Sd (w) 
can be found in a similar manner in the opposite limit d Sd , /d Sd and we have 




2w - b Sd - b s 



if a Sd /d Sd , -> (22a) 
if a Sd ,/a Sd ^0 (22b) 



Thus if the appropriate scale factor ratio vanishes for all d < d u , one can express (w) as a product of the parent 

distribution q and the i.i.d. functions T . Otherwise a novel distribution may be expected. If each integral in the 
product on the right hand side (RHS) of (|T6l) can be calculated, a further step is required for D > 1 to ascertain if 
there is a single scaling variable and the product is reducible to a single function. To illustrate these points, we now 
apply the above discussion to specific block fitness distributions: 

(i) For a class of block fitness distributions p(f) — je^f T /^ 1 " 7 , 7 > which decay as a power law, the location 
factor 1>n — and the scale factor — TV 1 / 7 so that the ratio a Sd /a Sd , = (sd/sd') 1 ^ 1 ■ As described below, two 
distinct cases arise depending on whether the number of mutations D is comparable to £: 

Case a. When D/£ — > as £ — > 00, using the approximation Sd ~ £ d /dl, we find that d Sd /d Sd , ~ g2d-D wn j cn 

vanishes for all d < d u for large I. Then due to (I22al) . the cumulative distribution V$ {w) for large £ can be written 
as 

"Kw-iWf)= n < 23 > 



d=Q 



where the cumulative Frechet distribution F^(y) — e y ,y > and zero otherwise [l|. In the limit w,£ — > 00 with 
W n — 2w/d Sn fixed, as the Frechet scaling function 

/ 2w \ _ ( W n a Sn \ { 1 if m < n (24a) 



a Sm ) \ a Sm J [ if m > n (24b) 

a nontrivial extreme value distribution is obtained if Wd = 2w/a SD is kept fixed. Thus we find that the maximum 
value distribution is given by 

^>)«fWf^) (25) 

which is verified in Fig. 2] for D = 5 and 7 = 1. As the above result holds for D = 1 as well (see (ITUl) ). we conclude 
that the universality class does not change from the i.i.d. class for 1 < D £ for parent distributions decaying as a 
power law. 
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FIG. 4: Scaled distribution for the first maximum for p(f) = f~ 2 e~ 1 ^ and D = 5 obtained using exact integration of (|16[1 for 
t = 40(+),60(x) and 80(o). The solid line shows the Frechet scaling function F w (y) = y- 2 e~ 1/y . 



Case b. For finite D/£ in the large £ limit, we can write 



Sd 
Sd> 



D-2d 

n 

n=l 



d + n 



exp 



r + x 



e-D + d + n 



1 - r 



R-lr 

£ I dx\n . 

o \l-R + r + x 

D-2d 

, R-2r = 



R-2r>0 



(26a) 
(26b) 



where r = and R = D/l. Since the integrand in the exponential on the RHS of (I26a[) is always negative, Sd/sd> 
decays to zero as £ increases for nonzero R — 2r but remains finite for vanishing R — 2r. As a result, the integral in 
(1171) or (|18l) does not reduce to i.i.d. distributions for d ^ d u and we may expect novel extreme value distributions 
for D ~ L 

(ii) A similar analysis can be carried out for bounded distributions p(f) = v{l — f) v_1 , v > for which 6jy = 1 
and a at = N^ 1 ^. When the number of mutations D — 1, a nonuniversal parent distribution is obtained by replacing 
£ by sd in (H2J) while a rescaled universal distribution holds for 1 < D <C I. Since the ratio d Sd ,/a Sd — > as f — > oo, 
due to (HI]) and ([22b]) . we get 



q(2w-l)l[P 



d=l 



(27) 



where the Weibull distribution IF^iy) — e ( ^ ,y < and unity otherwise Asa)->1 and i — >• oo, while the 
cumulative distribution g — >• 1, the distribution [2(w — l)/a Sm ] approaches zero for m > n and unity for m < n if 
M^n = 2(u> — l)/a s „ is kept fixed. Then a nontrivial maximum value distribution is obtained if 2(w — l)/a Sd remains 
finite and we obtain 



( 2{w-l) 



-F { 



(28) 



However for D ^ £, new extreme value distributions are expected for the reasons mentioned above. 

(iii) We finally consider block fitness distributions for which the scale factor ratio typically tends to a finite limit as 

£ — > oo. For p(f) = 5f 5 ~ 1 e~ i considered in Sec. IIII Al the scale factor ratio d Sd /d Sd , = (In Sd/ In Sd'Y 1 ^ 6 ^ 6 ,6^=1. 
Using Stirling's approximation for d > 0, we may write 



lnsrf rln((l — r)/r) — ln(l — r) 
lns d > ^ ? ln((l -r')/r') - ln(l - r') 



constant , if _D/£ — > 
, if d/l -+Q,D/£ finite 
constant , if both d/£,D/£ finite 



(29a) 
(29b) 
(29c) 



where r = d/£ and r' = d'/£ Thus we cannot express P^^w) in terms of the known functions even for 1 < D <C £ 
when 5^1. But some progress is possible for S — 1 which we discuss next. 
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As we have already mentioned, for p(f) = e~f the scale factor is independent of £ and thus the scale ratio is a 
constant for any D. However it is possible to find Q Sd for all d < d u , D > 1 for exponentially distributed fitnesses. 
Although the maximum value distribution for this case has been studied in a previous work [231 ] . here we present a 
simpler derivation. For large £, using (fTB")) for k = 1 in (fT5|) . we obtain 



2u< 



4fo(i - 

df e- fo 



e -^+fo) SDe -h e s D e-fo 



n 



2u< 



-2™ + / d _. _ 



see 



/o 



n 

d=l 



dfcte 



-fd e -s d s d ,e 



"ltd 



where we have used (|19[) for d = term. In the limits w, 
is finite and we get 



oo, the above integrals are nontrivial if Wd = s^s^e 



(30) 
(31) 

-2w 



dfoe- fa 1 - 



Ho 



W n 



V- 



n 

d=l 



dfde 



-fd e -W d /f d 



-W 



W o F(0, W )] 11 2^JW~ d K 1 (2^/w~ d ) 



(32) 
(33) 



d=l 



where ifi(a;) is the modified Bessel function of second kind [2J]. One can check that the above equation reduces to 
the cumulative distribution for D = 1 obtained in Sec. IIII Al on using Wo = e~" in (fT4|) [25]. We now consider the 
above cumulative distribution when D/t is zero and when it is finite in the large £ limit. 

Case a. Since the distribution Q Sd moves towards larger w as d increases (see Fig. [5]), we may try Wd u as a scaling 
variable to reduce the product in (|T5| to a single function. For D <C i on using « £ d /d! for large we see that 



W d d u \(D - d u 



d\{D-d)\ 



(34) 



which is finite for all d < d u . Thus unlike for the algebraically decaying or bounded block fitness distributions, all the 
factors in (|33|) contribute to the cumulative distribution V$ (w). 

Case b. For D ~ £, the ratio Wd/Wd u goes to zero for nonzero R — 2r and unity for vanishing R — 2r: 



W d 



n 



d + n 



D 



£-D + d + i 
I — d u +n 



exp 



/ r+x l—R+r+x 

dx In . — , ; x . . . 

\{R/2) + x l-(R/2)+x 



1 , R - 2r = 



, R - 2r (3Sa) 
(35b) 



Thus although each term in (|16l) can be calculated (unlike for other block fitness distributions discussed above) , most 
of the terms contribute to the distribution. As shown in Fig. [SJ the distribution Q Sd 's arc quite well separated for 
small d but get clustered at larger d which is consistent with the behavior of Wd/Wd u discussed above. 



IV. CONCLUSIONS 



A set of random variables obtained by summing a fixed number of i.i.d. random variables are correlated when they 
have at least one common term in the sum. Such nonindependent variables may describe diverse quantities such as 
breeding value of an animal [2(1 H3] , fitness of a protein or antibody [l7| and energy of a directed polymer [ll| . For 
two different linear combinations of i.i.d. random variables, we showed that as the sum is over a set of independently 
distributed variables, it is possible to write the extreme value distribution for weakly correlated random variables as 
an integral involving the maximum value distribution of i.i.d. random variables which makes the problem amenable to 
analysis at least in some cases. Interestingly even with weak correlations, a rich variety of extreme value distributions 
result: they can be highly nonuniversal parent distributions or universal extreme value distributions for i.i.d. random 
variables or novel distributions unrelated to these. 

When a single mutation occurs in the initial sequence, the limiting extreme value distribution (up to a rescaling) 
is found to be one of the i.i.d. distributions when the initial distribution p(/) decays faster than an exponential 
and parent distribution if p(f) decays slower than an exponential but a novel distribution appears for exponentially 
decaying p(f). Such a classification has also been observed in the context of near-extreme value statistics of i.i.d. 
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FIG. 5: Cumulative distribution (w) (bold line) for exponential parent distribution when D = 21 and (. = 80 obtained 
using (|30[) . The other curves show the individual factors Q Sd in the product in (|30[) with d — 0, d u (left to right). 



random variables |28( . The situation is more complex when multiple mutations are introduced in the initial sequence. 
When the number of mutations D does not scale with I so that D/l — »■ for large £, our analysis showed that D = 1 
and 1 < D C £ may exhibit different distributions. While the algebraically decaying parent distributions remain 
robust in that the universal Frechet distribution holds for 1 < D <C I, for bounded distributions the nonuniversal 
parent distribution for single mutation problem changed to universal Weibull distribution for D > 1. For unbounded, 
non-exponential parent distributions decaying faster than a power law also, the extreme value distribution for single 
mutation case fail to hold for D > 1. For all the three classes of block fitness distributions considered, novel extreme 
value distributions are expected if the number of mutations in both blocks is comparable to the block length. For 
exponentially distributed block fitnesses, we analysed the maximum value distribution for all D and gave explicit 
expressions for two new extreme value distributions (see (| 14[) and (|33[) 1. It would be very interesting to see if they 
occur in other extreme value problems as well. 

Acknowledgments: The author thanks S. Sabhapandit for pointing out reference [28j . 
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